N 


NASA Technical Memorandum 103100 


Nume ri cal Simulations of 
Supersonic Flow Through 
Oscillating Cascade Sections 


Dennis L. Huff 
Lewis Research Center 
Cleveland, Ohio 

and 


T.S.R. Reddy 

The University of Toledo 

Toledo, Ohio 


Prepared for the 

15th Southeastern Conference on Theoretical and Applied Mechanics 
sponsored by the Georgia Institute of Technology 
Atlanta, Georgia, March 22-23, 1990 


\ s t sm~ - ■ sw-.- - ^ ” •' Hfc • '■my**' 



IUAS/V 

( NAS A-TM-lOl 100) 
SUPERSONIC FLOW 

a sections (nasa) 


numerical simulations of 

[HROUGH OSCILLATING CAjCAOL 
10 r> CSCL 0 1 A 


NQ0-200SL 


G3/u? 


Unci js 
0775:346 



NUMERICAL SIMULATIONS OF SUPERSONIC FLOW 
THROUGH OSCILLATING CASCADE SECTIONS 


Dennis L. HufT 

National Aeronautics and Space Administration 
Lewis Research Center 
Cleveland, Ohio 44135 


and 

T.S.R. Reddy* 

The University of Toledo 
Toledo, Ohio 43606 


Abstract 

A finite difference code has been developed for modeling inviscid, unsteady supersonic flow by solution of the compressible 
Euler equations. The code uses a deforming grid technique to capture the motion of the airfoils and can model oscillating 
cascades with any arbitrary interblade phase angle. A flat plate cascade is analyzed, and results are compared with results 
from a small-perturbation theory. The results show very good agreement for both the unsteady pressure distributions and 
the integrated force predictions. The reason for using the numerical Euler code over a small-perturbation theory is the ability 
to model "real" airfoils that have thickness and camber. Sample predictions are presented for a cascade of loaded airfoils 
and show appreciable differences in the unsteady surface pressure distributions when compared with the flat plate results. 
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sonic velocity 
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blade chord length 

total energy of the fluid per unit volume 
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incidence angle 

Jacobian of transformation 

reduced frequency based on semichord, ~~ 
Mach number 

Reynolds number based on chord 
real part of { } 

arc length of a grid line in the redirection 
time normalized by — 

Cartesian velocities normalized by C i 
total velocity 

weighting function for grid deformation 


x , y Cartesian coordinates normalized by chord length 
a amplitude of pitching 

(3 flow angle 

Y stagger angle 

A CP pressure difference coefficient, CP _ - CP „ 
t\ normal direction of transformed 

coordinate system 

£ chordwise direction of transformed 

coordinate system 
p fluid density 

a interblade phase angle 

x time variable 

x ' thickness-to-chord ratio 

uo airfoil oscillation frequency 


Subscripts: 

1 conditions at the inlet 

2 conditions at the exit 

+ ,- upper and lower surfaces on airfoil, respectively 
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Introduction 

Many turbine engines operate with supersonic relative flow through the rotor and stator blade rows. Much research is 
needed to ensure stable operation of this type of design. Understanding the aeroelastic behavior and identifying the flutter 
boundaries is critical. Fundamental methods of flutter analysis depend on the successful predictions of blade loading and 
blade motion. Ideally, this is an interactive process where a structural analysis determines the blade motion from the blade 
loading and the loading is determined by the flow analysis from the blade motion. Many of the existing structural and 
aerodynamic analysis methods are based on flat plate theory and introduce approximations to model thickness and camber 
of the blades. One alternative is to use a finite difference algorithm to predict the flow field and account for these effects. 

The present research focuses on numerical solutions of two-dimensional supersonic flow through oscillating cascades. 
Numerical methods tend to converge faster in this flow regime compared with their subsonic counterparts. Harmonic pitching 
motions are prescribed for the blade sections for both zero and nonzero interblade phase angles. The code uses the deforming 
grid technique introduced in references 1 and 2 for convenient specification of the periodic boundary conditions from blade 
to blade. Several sample predictions are presented for oscillating cascades in this flow regime. The results and theory in this 
paper are based on the work presented in reference 3. This analysis can be applied to any fan or compressor designs with 
supersonic relative Mach numbers. 


Governing Equations 

A major portion of the present code is based on the unsteady solver developed by Sankar and Tang (ref. 4) for flow past 
isolated airfoils. This code solves the two-dimensional, unsteady, Reynolds-averaged, compressible Navier-Stokes equations 
in strong conservation form on a body-fitted moving coordinate system using an ADI procedure. These equations can be 
written as 






( 1 ) 


where 


(J * J~ ] {p , pu . pi/, e) 


( 2 ) 


and p is the fluid density, u and v are the Cartesian components of the fluid velocity, and e is the total energy per unit volume. 
The body-fitted ( £ , q f t) coordinate system is related to the Cartesian coordinates by using the following transformation: 

Tl-T](x,y t O (3) 

t - 1 

The Jacobian of the transformation is given by 


J = ^ T ly _T lx^ y 


1 

* t y n -* n y* 


and the metrics of the transformation are given by 


%y~~ Jx r\ fly"' 7 ** 


(4) 


(5) 


Standard central differences were used to compute x x , y t , x n , and y n which were then used to calculate the metrics. 

The P and £ terms in equation (1) are the inviscid terms in the % and redirections, respectively. The viscous terms, £and 
5, are omitted to give solutions of the Euler equations. In fact, all solutions presented in this paper are inviscid. The 
Beam-Warming block ADI algorithm is used to solve the governing equations. A Jameson-type artificial dissipation is added 
in both directions of the computational plane to increase solution stability. The solution is second-order accurate in space 
and first- or second-order accurate in time. Further information about the algorithm for isolated airfoils can be found in 
reference 4. 
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Grid 

A unique feature of the present code is the treatment of the grid for oscillating cascades. A method for deforming the 
grid was developed in references 1 and 2 for zero and nonzero interblade phase angles. The present study uses the same 
technique. 

The code uses a C-grid generated from Sorenson’s (ref.5) GRAPE code, which was modified by Chima (ref.6) for improved 
modeling in turbomachinery problems. One C-grid is generated for each blade in the cascade and undergoes a local 
deformation associated with the blade motion. The outer boundary of the C-grid is defined by the user in the GRAPE code. 
A deforming grid technique is used to locate the position of the grid as a function of time. The inner boundary moves with 
the blade, while the outer boundary remains fixed in space. The grid lines connecting the inner and outer boundaries are 
allowed to deform. The amount of deformation is a function of the distance away from the surface of the airfoil. A weighting 
function u;is defined as 






( 6 ) 


where s is the arc length of a grid line from the airfoil sur- 
face (r|= 1) to some grid point along £ = constant, and r| max 
is the outer boundary grid line. The grid deformation is defined as 

A^ (/ = ^(Ax ( ') 

a ® w iJ (Ay u ') (7) 

where A x 0 and Ay are the spatial differences that would exist between successive time steps if the entire grid were moved 
as a rigid body. From equations (6) and (7), we see that nodes at the inner boundary (s=0) give « 1, which means that 
the airfoil surface follows the rigid body motion of the blade. Conversely, the outer boundary nodes give w tJ * 0, and the 
node positions remain fixed at the initial specified locations. The interior nodes shear in space relative to the initial grid as 
io,-,- varies between 0 and 1. The node velocities can be easily found by dividing the grid deformation by the time step value. 

Multiple blade computations are made possible by stacking the C-grids for each blade and passing information between 
the periodic boundaries. Each C-grid is expanded by one grid line in the redirection at the outer boundary to provide ghost 
points for the implementation of the periodic boundary conditions. This allows the periodic boundary condition to be treated 
implicitly. The zero interblade phase angle is the simplest case for grid generation. Periodic boundary conditions are applied 
across the upper and lower boundaries and require a grid for only one blade. However, for supersonic axial flow, nonzero 
interblade phase angles require computations around two additional blades for exact treatment of the periodic boundary 
conditions. This increases the computational time by a factor of three, but provides an exact boundary condition. 

Figure 1 shows a deforming grid for a typical nonzero interblade phase angle computation. Multiblade solutions for 
oscillating cascades are done by generating grids for each blade and assembling them into a cascade. In this sample case, 
the interblade phase angle a is 180 degrees. A 199 by 22 C-grid is generated for an airfoil section. The solver automatically 
generates the grids for the adjacent blades before the unsteady solution begins by deforming the mean flow grid for one blade 
through one cycle of oscillation and saving the two grids that occur at the desired interblade phase angle. The grids are then 
assembled to give the initial multiblade grid to be used for the oscillating cascade solution at cot = 0 (fig. 1(a)). As time 
progresses in the unsteady solution, the grid for each blade deforms to model the specified interblade phase angle (fig. l(bl). 
Notice how the outer boundary of the grid around each blade remains fixed in space, while the inner boundary follows the 
motion of the airfoil. An actual run from an oscillating cascade typically has a small pitching amplitude and does not distort 
the grid as much as shown in figure 1. 


Boundary Conditions 

The present solution independently solves the flow equations for each grid around a blade and uses periodic boundary 
conditions along the upper and lower boundaries to model the cascade effects. Ghost points are assigned at the first interior 
gridline (q = q max - 1) and are used by the adjacent grid from the next blade. Although it is tempting to use the most current 
flow information as it becomes available from the integration scheme, it is important to only use flow information from the 
same time step across the periodic boundaries. This eliminates time inaccuracy due to the grid stacking direction of the 
multiple blade solutions. The metric data are also forced to be continuous along the periodic boundaries. This procedure 
essentially makes the blade-to-blade periodic boundaries invisible to the flow solution. 

In supersonic flow, the domain ot the flow field can be reduced by inspecting the shock structures through the cascade. 
The disturbances generated inside the Mach cone of the leading edge do not influence the flow outside the cone, which means 
that only three blades need to be considered in the solutions. For example, consider the cascade geometry defined in figure 
2, and let the flow conditions be defined such that the bow shock off the leading edge intersects the adjacent blade surfaces. 
The flow above blade 2 and below blade 3 will have no influence on the surface pressures on blade 1. For this reason, the 
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boundary conditions on the upper and lower boundaries of the global grid (the grid containing all three C-grids) can be 
arbitrary. The present analysis allows the upper and lower boundaries of the global grid to be periodic, which means that 
only the surface pressures from the reference blade will be valid for the specified interblade phase angle. 

The inlet conditions are assumed to be uniform by specifying the flow density, velocity, flow angle, and pressure. The exit 
flow variables are extrapolated from the interior by using a simple first-order model. These boundary conditions are valid 
for supersonic flow at the inlet and exit. Solid wall boundary conditions are applied along the airfoil surface, and the flow 
variables are averaged across the slit aft of the airfoil. 


Results and Discussion 

Sample predictions for an oscillating cascade in supersonic flow have been done for both a cascade of flat plates and a 
cascade of airfoils suitable for supersonic flow. There are no unsteady experimental data available for comparisons with the 
airfoil section predictions. At present, any validation of this code is limited to comparisons with a small-perturbation theory 
for flat plates. 


Flat Plate Cascade 

Before an unsteady solution can be obtained for an oscillating cascade, a good steady-state solution for the mean flow 
conditions should be calculated. As previously described, a C-grid is generated for one blade by using a modified GRAPE 
code. For the present solutions, a 199 by 22 grid in the and r> directions, respectively, is used around a flat plate with the 
thickness-to-chord ratio, x =0.005, The leading and trailing edges were rounded to aid in C-grid generation, and therefore 
this is only an approximate representation of a flat plate. Three C-grids are then assembled to define the geometry for the 
cascade, as shown in figure 3. The present test case was selected from Kielb and Ramsey (ref, 7), who simulated the rotor 
design for the supersonic throughflow fan by using Lane’s theory (ref, 8) to determine the unsteady pressure distributions. 
The stagger angle y is 28 degrees and the gap-to-chord ratio g / cis 0.311. An inviscid, steady-state solution has been generated 
using the present code for M ^2.61 and (3 j = 2 8 degrees, which gives an incidence angle i of 0 degrees. The steady-state 
solution required about 1500 iterations with Ax=0.005 for the maximum residual to drop four orders of magnitude, which 
corresponds to about 100 sec of CPU time on the CRAY-YMP. This time can be reduced by implementing a variable local 
time step instead of the constant time step used in the present solution. Since the focus of the present research is the solution 
of the unsteady Euler equations, little time has been spent on accelerating the steady-state solutions. 

Solutions were done for oscillating flat plates for various interblade phase angles and were compared with the small- 
perturbation theory of Lane (ref. 8). The unsteady solutions were run with first-order temporal accuracy. The surface pressure 
time histories were recordea and found to reach a reasonably periodic solution after three cycles of cascade oscillation. A 
Fourier transform was done on the fourth cycle to determine the first harmonic pressure distribution relative to the airfoil 
motion. The pressures were normalized by the airfoil pitching amplitude and the phase was referenced to the airfoil pitching 
motion starting at the maximum (nose-up) blade angle. Predictions for the unsteady pressure difference distributions (A C P) 
on the flat plates are shown in figure 4 for o=0 and 180 degrees, fc=0.50, and a=0.10 degrees. The small amplitude of 
oscillation is used for comparisons with the small-perturbation theory. The agreement with Lane’s theory is good for both 
the real and imaginary parts of pressure, except that the locations of the shock reflections were predicted about 10 percent 
forward of the Lane’s predictions. The actual locations of the shocks are exact for an infinitely thin flat plate, as used by 
Lane’s theory. The finite thickness approximation of the fiat plate used in the Euler code is expected to be the greatest 
contributor to this error. The wiggles near the leading edge of the Euler solution are numerical and are reduced for real 
airfoils. It is possible to adjust the numerical dissipation coefficients to minimize these wiggles; however the shocks tend to 
smear, and the accuracy of the solution becomes jeopardized. The force coefficients obtained by integrating these pressures 
are insensitive to these wiggles. 

The unsteady solutions use 2.11 x 1 0 " 5 sec of CPU per time step per grid point per blade and required 2047 time steps to 
complete 4.25 cycles of oscillation. In practice, only three cycles of oscillation are necessary to reach a periodic solution, 
which corresponds to 400 sec of CRAY-YMP CPU time. 


Airfoil Section Cascade ,_ r „ 

The reason for using the Euler code over Lane’s theory is the ability to model "real" blade sections that have thickness 
and camber. While Lane’s theory offers the most practical analysis tool for investigating a wide range of flutter parameters, 
the Euler code can be used to check the effects of blade loading. The cascade of airfoils shown in figure l jijas been selected 
to investigate the real blade effects. The parameters for the flow and cascade geometry are similar to those used for the flat 
plate analysis, except that the airfoils are now loaded: g / c= 0.311, degrees, M i = 2.61, and p i = 36.0 degrees. The grid 
presented in figure 1 is used in the present investigation and has the same dimensions and parameters described for the grid 
used in the flat plate analysis. 
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Two unsteady solutions for the blades oscillating with o= 0 degrees and 180 degrees, £=0.50, and a=0.1 degrees are 
presented in figure 5. The discontinuities near 20- and 70-percent chord are evidence of the shock structures. The shock 
movement due to a pitching amplitude of 0.1 degrees was expectedly small. The addition of thickness and camber to the 
solution caused the shocks to move forward from their locations in the flat plate analysis. 

Conclusions 

A code has been developed that solves the nonlinear flow field for oscillating cascades in supersonic flow and can be used 
in flutter analysis. A finite difference code has been developed for modeling compressible, inviscid, unsteady supersonic flow 
by solution or the Euler equations. The code uses a deforming grid technique to capture the motion of the airfoils and can 
model oscillating cascades with any arbitrary interblade phase angle. A fiat plate analysis is done for comparisons with a 
small-perturbation theory. The results show good agreement for the unsteady pressure distributions predictions. Sample 
predictions are presented for a loaded cascade and demonstrate how different the unsteady pressure distributions can be 
when including the effects of thickness and camber. 

The code is a tool for modeling '’real" blade sections, which can have significantly different flow characteristics than results 
from methods that are restricted to flat plate geometries. Obviously, further research is needed to validate the code, such as 
detailed experimental data for oscillating cascades. It is now possible to add a simple two-dimensional structural model to 
evaluate the effects of structural properties. 
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(a) ut = 0°. 


(b) ut = 90°. 

Figure t .—Deforming grid technique with exaggerated motion. 
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(a) Global grid. 

(b) Enlarged view of leading-edge region. 
Figure 3.— Grid for flat plate cascade. 
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Figure 4 ,— First harmonic pressure difference distribution for flat plate cascade. Stagger angle, 7 , 28 # ; gap-to-chord ratio, g/c, 0.311; inlet flow angle, 
0,, 28 # ±0.10 a ; reduced frequency based on semichord, k, 0.5. 





RATIO OF NORMALIZED LENGTH TO CHORD, X'/C 


Figure 5.— First harmonic pressure difference distribution for rotor cascade. Inlet Mach number, M„ 2.61; stagger angle, y, 28“; gap-to-chord ratio, 
g/c, 0.311; inlet flow angle, 36“=fc0.10“; reduced frequency based on semichord, k, 0.5. 
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